#include <fstream>
#include <iostream>
#include <iomanip>
#include "TFile.h"
#include "TTree.h"
#include "TSystem.h"
#include "TH1.h"
#include "TH2.h"
#include "TEllipse.h"
#include "TCanvas.h"
#include "TGraph.h"
#include "TStyle.h"
#include "TLegend.h"
#include "TClonesArray.h"

using namespace std;

void Neutron(){

  int ret = gSystem->Load( "../Objects/libObjectClass.so" );
  if( ret!=0 ){
    cout << "Failed to Load Library: ret = " << ret << endl;
    return;
  }

  TFile *fread=new TFile("../outputFiles/test.root");
  //  TFile *fread=new TFile(Form("../outputFiles/Ar_0V_%dkeV.root",ene));
  TTree *tree=(TTree*)fread->Get("tree");
  
  TClonesArray *JPMTa=new TClonesArray("PMT");
  TClonesArray *MPMT1a=new TClonesArray("PMT");
  TClonesArray *MPMT2a=new TClonesArray("PMT");
  TClonesArray *MPMT3a=new TClonesArray("PMT");
  TClonesArray *MPMT4a=new TClonesArray("PMT");
  int nPhoton;

  tree->GetBranch("JPMTb")->SetAutoDelete(kFALSE);
  tree->GetBranch("MPMT1b")->SetAutoDelete(kFALSE);
  tree->GetBranch("MPMT2b")->SetAutoDelete(kFALSE);
  tree->GetBranch("MPMT3b")->SetAutoDelete(kFALSE);
  tree->GetBranch("MPMT4b")->SetAutoDelete(kFALSE);
  tree->SetBranchAddress("JPMTb",&JPMTa);
  tree->SetBranchAddress("MPMT1b",&MPMT1a);
  tree->SetBranchAddress("MPMT2b",&MPMT2a);
  tree->SetBranchAddress("MPMT3b",&MPMT3a);
  tree->SetBranchAddress("MPMT4b",&MPMT4a);
  tree->SetBranchAddress("nPhoton",&nPhoton);

  TH1F *hNphoton=new TH1F("hNphoton","hNphoton",5000,0,5000);
  TH1F *hNJhit[18];
  TH1F *hNM1hit[18];
  TH1F *hNM2hit[18];
  TH1F *hNM3hit[18];
  TH1F *hNM4hit[18];

  TLegend *leg=new TLegend(70,170,95,240,"","user");

  for(int i=0;i<18;i++){
    hNJhit[i]=new TH1F(Form("J_%d",i),Form("J_%d",i),100,0,100);
    hNM1hit[i]=new TH1F(Form("M1_%d",i),Form("M1_%d",i),100,0,100);
    hNM2hit[i]=new TH1F(Form("M2_%d",i),Form("M2_%d",i),100,0,100);
    hNM3hit[i]=new TH1F(Form("M3_%d",i),Form("M3_%d",i),100,0,100);
    hNM4hit[i]=new TH1F(Form("M4_%d",i),Form("M4_%d",i),100,0,100);

    hNM1hit[i]->SetLineColor(2);
    hNM2hit[i]->SetLineColor(3);
    hNM3hit[i]->SetLineColor(4);
    hNM4hit[i]->SetLineColor(6);

    if(i==0){
      leg->AddEntry(hNJhit[i],"J-PMT","l");
      leg->AddEntry(hNM1hit[i],"Meg1","l");
      leg->AddEntry(hNM2hit[i],"Meg2","l");
      leg->AddEntry(hNM3hit[i],"Meg3","l");
      leg->AddEntry(hNM4hit[i],"Meg4","l");
    }

  }


  int nev = tree->GetEntries();
  cout<<"NEV = "<< nev <<endl;

  int counter = 0;
  int counter1 = 0 ;
  int counter2 = 0 ;
  int counter3 = 0 ;
  int counter4 = 0 ;

  for(int iev=0;iev<nev;iev++){
    JPMTa->Clear();
    MPMT1a->Clear();
    MPMT2a->Clear();
    MPMT3a->Clear();
    MPMT4a->Clear();
    tree->GetEntry(iev);

    int iz = iev / 10000;
    float zpos = -545.+iz*10.;

    hNphoton->Fill(nPhoton);

    int njhit = JPMTa->GetEntriesFast();
    int nmhit1= MPMT1a->GetEntriesFast();
    int nmhit2= MPMT2a->GetEntriesFast();
    int nmhit3= MPMT3a->GetEntriesFast();
    int nmhit4= MPMT3a->GetEntriesFast();

    if( nPhoton>0 ){
      hNJhit[iz]->Fill( njhit );    
      hNM1hit[iz]->Fill( nmhit1 );    
      hNM2hit[iz]->Fill( nmhit2 );    
      hNM3hit[iz]->Fill( nmhit3 );    
      hNM4hit[iz]->Fill( nmhit4 );    

      counter++;
      if(iz==0) counter1++;
      if(iz==6) counter2++;
      if(iz==12) counter3++;
      if(iz==17) counter4++;
    }

  }

  int Jint[18]={0};
  int M1int[18]={0};
  int M2int[18]={0};
  int M3int[18]={0};
  int M4int[18]={0};
  for(int iz=0;iz<18;iz++){
    for(int ibin=51;ibin<102;ibin++){
      Jint[iz] += hNJhit[iz]->GetBinContent(ibin);   
      M1int[iz]+= hNM1hit[iz]->GetBinContent(ibin);
      M2int[iz]+= hNM2hit[iz]->GetBinContent(ibin);
      M3int[iz]+= hNM3hit[iz]->GetBinContent(ibin);
      M4int[iz]+= hNM4hit[iz]->GetBinContent(ibin);
      if(iz==0&&(ibin==101||ibin==100)) cout<<hNJhit[iz]->GetBinContent(ibin)<<endl;
    }
  }

  TGraph *gr1=new TGraph();
  TGraph *gr2=new TGraph();
  TGraph *gr3=new TGraph();
  TGraph *gr4=new TGraph();
  TGraph *gr5=new TGraph();
  for(int i=0;i<18;i++){
    gr1->SetPoint(i,i+2.5,(double)Jint[i]);
    gr2->SetPoint(i,i+2.5,(double)M1int[i]);
    gr3->SetPoint(i,i+2.5,(double)M2int[i]);
    gr4->SetPoint(i,i+2.5,(double)M3int[i]);
    gr5->SetPoint(i,i+2.5,(double)M4int[i]);
  }
  TCanvas *cc=new TCanvas("c","c",800,600);
  cc->SetGrid(1);
  TH1 *frame=gPad->DrawFrame(0,0,20,2000);
  frame->SetXTitle("Distance from Cathode (cm)");
  frame->SetXTitle("# of Photon");
  frame->Draw();
  gr1->SetLineColor(1);
  gr2->SetLineColor(2);
  gr3->SetLineColor(3);
  gr4->SetLineColor(4);
  gr5->SetLineColor(6);
  gr1->SetMarkerStyle(21);
  gr2->SetMarkerStyle(21);
  gr3->SetMarkerStyle(21);
  gr4->SetMarkerStyle(21);
  gr5->SetMarkerStyle(21);
  gr1->Draw("plsame");
  gr2->Draw("plsame");
  gr3->Draw("plsame");
  gr4->Draw("plsame");
  gr5->Draw("plsame");

  cout<<"#(1>photon)          = "<<counter<<endl;
  cout<<"#(&& Cathode+ 2.5cm) = "<<counter1<<endl;
  cout<<"#(&& Cathode+ 8.5cm) = "<<counter2<<endl;
  cout<<"#(&& Cathode+14.5cm) = "<<counter3<<endl;
  cout<<"#(&& Cathode+19.5cm) = "<<counter4<<endl;


  gStyle->SetOptStat(000000);

  TCanvas *c1=new TCanvas("c1","c1",800,600);
  c1->SetLogy();
  hNphoton->SetXTitle("# of Photon");
  hNphoton->SetYTitle("NEV");
  hNphoton->Draw();

  TCanvas *c2=new TCanvas("c2","c2",1200,1000);
  c2->Divide(2,2);
  c2->cd(1);
  TH1 *hframe=gPad->DrawFrame(0,0,100,250);
  hframe->SetXTitle("# of Photon");
  hframe->SetYTitle("NEV");
  hframe->SetTitle(Form("Cathode + %.1fcm",2.5));
  hframe->Draw();
  hNJhit[0]->Draw("same");
  hNM1hit[0]->Draw("same");
  hNM2hit[0]->Draw("same");
  hNM3hit[0]->Draw("same");
  hNM4hit[0]->Draw("same");
  leg->Draw();
  c2->cd(2);
  TH1 *hframe1 = (TH1*)hframe->Clone();
  hframe1->SetTitle(Form("Cathode + %.1fcm",8.5));
  hframe1->Draw();
  hNJhit[6]->Draw("same");
  hNM1hit[6]->Draw("same");
  hNM2hit[6]->Draw("same");
  hNM3hit[6]->Draw("same");
  hNM4hit[6]->Draw("same");
  leg->Draw();
  c2->cd(3);
  TH1 *hframe2 = (TH1*)hframe->Clone();
  hframe2->SetTitle(Form("Cathode + %.1fcm",14.5));
  hframe2->Draw();
  hNJhit[12]->Draw("same");
  hNM1hit[12]->Draw("same");
  hNM2hit[12]->Draw("same");
  hNM3hit[12]->Draw("same");
  hNM4hit[12]->Draw("same");
  leg->Draw();
  c2->cd(4);
  TH1 *hframe3 = (TH1*)hframe->Clone();
  hframe3->SetTitle(Form("Cathode + %.1fcm",19.5));
  hframe3->Draw();
  hNJhit[17]->Draw("same");
  hNM1hit[17]->Draw("same");
  hNM2hit[17]->Draw("same");
  hNM3hit[17]->Draw("same");
  hNM4hit[17]->Draw("same");
  leg->Draw();



}
